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In this work it is shown that the use of a hydrodynamical model of heavy ion collisions which 
incorporates recent developments, together with updated photon emission rates, greatly improves 
agreement with both ALICE and PHENIX measurements of direct photons, supporting the idea 
that thermal photons are the dominant source of direct photon momentum anisotropy. The event- 
by-event hydrodynamical model uses IP-Glasma initial states and includes, for the first time, both 
shear and bulk viscosities, along with second order couplings between the two viscosities. The effect 
of both shear and bulk viscosities on the photon rates is studied, and those transport coefficients 
are shown to have measurable consequences on the photon momentum anisotropy. 


I. INTRODUCTION 


The collision of heavy ions is the only way to produce 
and study hot and dense strongly interacting matter in 
the laboratory. Therefore, relativistic nuclear collisions 
represent a unique opportunity to explore QGD (Quan¬ 
tum ChromoDynamics) - the theory of the strong inter¬ 
action - in extreme conditions of temperature and den¬ 
sity. A vibrant experimental program is currently under 
way at RHIG (the Relativistic Heavy Ion Collider, at 
Brookhaven National Laboratory, NY) and at the LHC 
(the Large Hadron Collider, at CERN, in Geneva). The 
data accumulated at these facilities represent unequivo¬ 
cal evidence that a new state of matter has been created 
in the collision of large nuclei: the Quark-Gluon Plasma 
(QGP) [1]. One of the most striking features of the QGP 
is that its time-evolution can be modelled with relativis¬ 
tic fluid dynamics. Early studies concentrated on ideal 
fluids [2] , but the realization that hadronic data from rel¬ 
ativistic heavy ion collisions could be used to extract the 
transport coefficients of QGD - in particular the shear 
viscosity to entropy density ratio, ??/s - opened new hori¬ 
zons. More specifically, the azimuthal anisotropy of the 
particle momentum distribution in a single event can be 
characterized by Vn, the coefficients of Eourier expansion 
in azimuthal angle 4> [3]: 
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where px is the transverse momentum, the A/niPT, y) are 
orientation angles, and y is the rapidity. Viscous hy¬ 
drodynamics calculations have established a quantitative 
connection between the empirically-extracted u„’s and 
the value of 77 /s, the shear viscosity to entropy density 
dimensionless ratio [4, 5]. Studies of the different exper¬ 
imental anisotropic flow coefficients have concluded that 


the phenomenologically-extracted rj/s values were close 
to the conjectured lower bound of 77 /s = l/ 47 r [ 6 ]. The 
ability for fluid-dynamical models to quantitatively re¬ 
produce the measured behavior of the flow anisotropy 
coefficients of hadrons has been one of the major high¬ 
lights of the entire relativistic heavy-ion program [7]. 

A precise determination of quantities such as 77 /s is 
made difficult by significant uncertainties in the descrip¬ 
tion of the early time dynamics and the effect of addi¬ 
tional sources of dissipation [8], among others. Electro¬ 
magnetic observables produced in ultrarelativistic heavy 
ion collisions can be used as additional probes of the 
properties of the QGP, and can help constrain the trans¬ 
port coefficients of QGD. The task of measuring photons 
and subtracting the large background of hadronic decay 
photons has been undertaken at RHIG [9-12] and the 
LHG [13-15], and the direct photon transverse momen¬ 
tum spectra and azimuthal anisotropy are available at 
both colliders. The observation that the magnitude of the 
direct photon V 2 was similar in size to that of hadrons, 
along with the exponential behaviour of the measured 
direct photon spectra at low transverse momentum, sug¬ 
gested that both photons and hadrons were produced by 
a similar mechanism: (quasi-)thermal production. 

While event-by-event hydrodynamical models of heavy 
ion collisions were shown repeatedly to provide a good 
description of hadronic observables [16, 17], similar at¬ 
tempts at describing direct photon measurements did 
not meet with the same success [18]. A simultaneous 
description of the direct photon spectra and momentum 
anisotropy proved to be a particular challenge. In re¬ 
sponse to this apparent tension with measurements, in¬ 
vestigations into additional photon production mecha¬ 
nisms multiplied (e.g. [19-25]). 

In this paper, a hydrodynamical calculation of direct 
photon production is presented. It uses an up-to-date 
hydrodynamical model of heavy ion collisions [8] along 
with the latest photon emission rates [26-28]. Empha- 
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FIG. 1. The temperature dependence of the bulk viscosity 
to entropy density ratio, as used in this work. The points 
on the low temperature side are the results of a calculation 
from Ref. [32], while the high-temperature results are from 
Ref. [33]. 


sis is put on photon emission from the expanding QCD 
medium, referred to as “thermal photons”. 

The aim of this work is to present an up-to-date calcu¬ 
lation of thermal photons using the latest developments 
in hydrodynamical simulation of heavy ion collisions. Un¬ 
derstanding the current status of thermal photon produc¬ 
tion in heavy ion collisions will help guide future effort at 
identifying and constraining alternative photon produc¬ 
tion mechanisms. 


II. HYDRODYNAMICAL MODEL 

The relativistic fluid dynamics background that pro¬ 
vides the time-dependent environment in which the 
photon-generation mechanisms evolve, is the same here 
as that used for hadrons in Ref. [8]. We summarize its 
main features again here, for convenience. The initial 
state of the nuclear collision are modelled using the IP- 
Glasma approach [29], which builds on the “impact pa¬ 
rameter dependent saturation model” (IP-Sat) [30, 31] 
that constrains the distribution of initial colour sources 
drawing from electron-proton and electron-nucleus colli¬ 
sion data. The gluon fields are then evolved in space and 
time using classical Yang-Mills equations: [H^, F^"^] = 0, 
up to a proper time tq of order of the inverse of the sat¬ 
uration scale. The energy density, e, and the flow veloc¬ 
ities, u^, from the Yang-Mills evolution are then used to 
initialize the hydrodynamical evolution. This is achieved 
by solving u^(ro)r^YM('^o) = i{ro)u''{To) where is 

the classical Yang-Mills energy-momentum tensor. As in 
Ref. [8], To = 0.4 fm is used in this work. The IP-Glasma 
initial conditions are boost-invariant, and the subsequent 
hydrodynamical evolution is 2 -|- ID as well. 

The hydrodynamic evolution involves a dissipative 


term in the stress-energy tensor: 

, ( 2 ) 

where — u^u’'. In the above, = 

diag{l,—l,—l,—l) is the Minkowski tensor, is the 
shear-stress tensor, and 11 is the bulk pressure term. The 
time-evolution of these last two quantities is obtained by 
solving relaxation-type equations [34, 35]: 

Tnli + n = —-I-, (3) 
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is the double, symmetric, and traceless projection oper¬ 
ator. The expansion rate of the fluid \s 6 = and 

the shear tensor . In the present work, the 

shear-stress tensor and the bulk pressure are initialized 
to zero at time tq. 

The second-order transport coefficients rn, ^nn, Autt, 
t-k, Vi ^ 7 , Trir, and A^n are related to the shear vis¬ 
cosity g and bulk viscosity Q using formulae derived from 
the Boltzmann equation near the conformal limit [35]. 

Importantly, the hydrodynamic evolution stage is fol¬ 
lowed by a phase where discrete particles are pro¬ 
duced through the Cooper-Frye procedure [36]. Late 
stage hadrons further interact and freeze-out dynamically 
through the UrQMD approach and algorithms [37]. 

Since viscous hydrodynamics is used, the medium is 
not exactly in thermal equilibrium. Consequently, when¬ 
ever specific particle distributions are invoked - in the 
Cooper-Frye scheme or for thermal photon production 
- these will receive viscous corrections. In the present 
work, the momentum distribution fg/piPj X) is derived 
in Appendices (A) and (B), and is given by 

fs/FiP, X) = 4%(P) + X) + 5/B“/^(P, A) , 

( 5 ) 
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with as = 1 for bosons and ap = ~1 for fermions, with 
f^pjp^P) being correspondingly either the Fermi-Dirac or 
Bose-Einstein distribution. The pressure V, energy den¬ 
sity e, flow velocity rt, speed-of-sound and temperature 
T entering into the distribution functions are evaluated at 
spacetime point X. The momentum of the quasi-particle 
of mass m is denoted P. The bulk relaxation time is rn- 

Within the hybrid approach used here (IP-Glasma - 
dissipative hydrodynamics - UrQMD) a recent analysis 
of ALICE and CMS measurements indicates that identi¬ 
fied particle spectra, multiplicity distributions, and mul¬ 
tiple flow coefficients (u„{2}, n = 2,3, 4) can be globally 
reproduced, for pions, kaons, and protons [8]. The cal¬ 
culation and data analyses were done for 0 - 5% through 
30 - 40% centrality classes. These LHC-energy analy¬ 
ses, which include both bulk and shear viscosity, lead 
to rj/s = 0.095. A relatively narrow C,/s temperature- 
profile, illustrated in Fig. 1, was used. This profile peaks 
such that C/'S(Tpeak) ~ 0.3, with Tpeak = 180 MeV. As 
a reminder, this is a somewhat novel feature that most 
fluid-dynamical approaches to the modelling of relativis¬ 
tic heavy-ion collisions do not yet include. The approach 
used here also notably features non-linear terms that cou¬ 
ple the shear and bulk sectors of the viscous hydrodynam¬ 
ics [35]. 

Describing the late stage dynamic of the medium with 
an afterburner leads to significant complications in the 
evaluation of photon emission. As a consequence, late 
photon emission is not evaluated with the afterburner, 
but rather with hydrodynamics. The exact approach 
used is explained later in this work. 
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III. PHOTON SOURCES 


The photons measured in relativistic nucleus-nucleus 
collisions come from a variety of different sources, and 
these will be discussed in turn in this section. They fall 
in two broad categories: those with a thermal origin and 
those coming from “cold” processes. In searches for sig¬ 
nals from the quark-gluon plasma, the latter are usually 
thought of as a background. Since this background is ir¬ 
reducible in the experimental measurements, they never¬ 
theless deserve our full attention. We start by describing 
how prompt photon production is evaluated in this work. 


A. Prompt photons 

In the very first instants of the nuclear collision, the 
interacting nucleons will produce photons through par- 
tonic Compton interactions and quark-antiquark anni¬ 
hilations. In addition, QCD jets will be generated and 
these jets will fragment into many final states, some of 
which will include photons. 

The calculation of photon production in hadronic in¬ 
teractions using the techniques of perturbative QCD has 


FIG. 2. Top panel: Direct photon spectrum measured in 
^ysNN = 200 GeV proton-proton collisions at RHIG com¬ 
pared with perturbative QGD calculations at different scales 
Q. Bottom panel: normalised perturbative QGD calculations; 
see details in the main text. 


a long history [38] which has led to a fairly mature un¬ 
derstanding of the subject. The photon production cross 
section in proton-proton collisions can be written con¬ 
cisely as 


(Pp 


^ ^ fa/p (^a: Qfact) ^ fb/p Qfact) 

a,b,c,d 

dd" {Qren) ^ ^7/c Qfrag) 5 (8) 


where Qfact, Qren, and Quag are energy scales entering 
respectively into the parton distribution function fa, the 
partonic cross-section da, and the fragmentation func¬ 
tion D^/c- The cross-section da{Qren) is evaluated as a 
perturbative expansion in the strong coupling constant 
as{Q), and Qren is the scale at which as{Q) is evalu¬ 
ated. This scale, along with the factorisation and frag¬ 
mentation scales, should typically be of the order of the 
transverse momentum of final state partons. 

Computing photon production in proton-proton 
collisions from the formalism described above 
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RHIC Au-Au ^snn = 200 GeV 

LHC Pb-Pb ^snn = 2760 GeV 

Gentrality 

0-20% 

20-40% 

0-20% 

20-40% 

(Neon) 

793 

323 

1231 

501 


TABLE I. Average number of binary collisions in different centrality classes, at RHIC and the LHC. Nucleon positions sampled 
from a Wood-Saxon distribution, which are input of the IP-Glasma model, are used to evaluate the number of binary collisions 
in each event with the MC-Glauber model. Centralities are defined using the gluon multiplicity of each event, as described in 
Ref [39]. 


requires a proton parton distribution function 
fa/p i^a^Q fact) i a parton-to-photon fragmentation 
function D^/f,(zc,Qhag), and the partonic cross-section 
d(T((5ren)- The latter is currently known at next-to- 
leading order in the strong coupling constant for both 
isolated photons [40] and fragmentation photons [41]. 
Combined with next-to-leading order parton distribution 
and fragmentation functions, photon production using 
perturbative QCD has been shown to agree very well 
with direct photon measurements in proton-proton 
collisions at RHIC, at the LHC, and at previous 
colliders [42]. 

At high pj, prompt photons are by far the dominant 
source of direct photons. Those are calculated in heavy 
ion collisions by multiplying the number of photons pro¬ 
duced in proton-collisions by the number of binary col¬ 
lisions [43, 44]. The scaling procedure may be applied 
either to a ht of direct photon measurements in proton- 
proton collisions, or to a perturbative QCD calculation 
of photon production. The latter is used in this work, in 
conjunction with nuclear parton distribution functions 
EPS09 [45], which take into account cold nuclear matter 
effects. The study of fragmentation photon energy loss 
and jet-medium photon productions, two effects that are 
understood to modify low pj. prompt photon production 
in heavy ion collisions, is not undertaken here and will 
be the subject of a separate work. 

Prompt photons are added to other sources of direct 
photons on an event-by-event basis. The perturbative 
QCD calculation is thus scaled by the number of binary 
collisions in each event individually. For reference, we 
quote in Table I the centrality-averaged number of binary 
collisions in each centrality studied in this work. 

The pQCD framework used in this work is essentially 
the next-to-leading order calculation contained in the nu¬ 
merical code INCNLO [46]. The proton parton distribu¬ 
tion function and photon fragmentation function used are 
respectively CTEQGlm [47] and BFG-2 [48]. The factori¬ 
sation, renormalisation and fragmentation scales are all 
set equal to each other: Qfact = Qren = Qfrag = Q- The 
transverse momentum of the produced photon is used to 
set the scale, with a normalization constant Q = 

The effect of changing the proportionality constant be¬ 
tween Q and pf is essentially a change in the normal¬ 
ization of the prompt photon spectra. This can be seen 
in Fig. 2. The top panel shows the perturbative calcu¬ 
lation of prompt photons in proton-proton collisions for 
different choices of scale Q, from Q = pfp/2 to Q = 8p^. 


The lower panel shows the same calculations scaled by 
a constant so that they have the same normalization. It 
is clear from this last figure that the calculations overlap 
very well, showing that they have the same p^ depen¬ 
dence. It was verified that changing the scale Q at LHC 
energies also has the same effect on the photon spectrum, 
i.e. it changes the normalisation but not the transverse 
momentum dependence of the calculation. 

It is apparent that a small proportionality constant 
between Q and p'f, such as Q = pJ/2, provides a better 
description of the available measurements at RHIC. On 
the other hand, it can be seen in Fig. 2 that calculations 
are limited to > (1.5 CeV)/A, where A is the pro¬ 
portionality constant between Q and p'f. This limitation 
results from the presence of a scale Qo ~ 1.5 CeV, which 
is typically taken as the limit of applicability of perturba¬ 
tive QCD. Parton distribution functions and fragmenta¬ 
tion functions are usually limited to Q > Qo ~ 1-5 CeV. 
Calculations made with INCNLO are also subject to this 
limit in Q. Although this might appear to limit the value 
of py at which prompt photons can be evaluated with 
perturbative QCD, the scaling behaviour observed on the 
bottom panel of Fig. 2 shows that it is not the case: the 
effect of the scale Q is simply a change in normalization, 
and prompt photons can be evaluated at low pj by us¬ 
ing e.g. Q = 4py and changing the normalization of the 
calculation to that of Q = pJ/2. This is the procedure 
adopted in this work. 

It remains that pQCD is based on the idea that a large 
momentum exchange occurs in a hadronic collisions, al¬ 
lowing for a part of the cross-section to be computed 
perturbatively. It is understood that perturbative QCD 
eventually breaks down at low transverse momentum, al¬ 
though the exact value of pj. at which this happens is not 
clear. As discussed above, the perturbative QCD calcu¬ 
lation of prompt photons is in good agreement with the 
low py direct photon measurements in proton-proton col¬ 
lisions, although it slightly overestimates the very lowest 
point around 1 CeV. It remains to be seen how much 
lower in transverse momentum the agreement with data 
persists. 

While there are currently no low p^ photon data from 
proton-proton collisions at the LHC, the same extrapo¬ 
lation procedure was used to evaluate low px tt® produc¬ 
tion as an additional verification of the approach. Cood 
agreement with measurements was again found down to 
PT -- 1 - 2 CeV [49]. 

It is worth noting that, besides clarifying the domain of 
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avoids possible double-counting issues between mesonic 
and baryonic contributions. Finally, this work includes 
also recent estimates of tttt bremsstrahlung contributions 
[27], and of the reactions np —>■ wy, ttuj —>■ py, and 
TTUJ —>■ pn [28], absent from Ref. [26]. It is instructive 
to compare rates, prior to integrating them with a dy¬ 
namical four-volume evolution. This is done in Fig. 3. 
The figure shows the LO partonic rates of Ref. [52] (solid 
lines) compared with the hadronic rates of Refs. [26-28] 
(dashed lines) for a range of temperatures in the cross¬ 
over region. 


C. Non-cocktail hadronic decay photons 


FIG. 3. Ideal QGP and hadronic photon rate near the cross¬ 
over region. 

validity of perturbative QCD calculation of prompt pho¬ 
tons, additional measurements will also help constrain 
uncertainties due to the photon fragmentation function, 
which are significant in the soft domains of perturbative 
QCD calculation of prompt photons [50]. Without direct 
measurements, those uncertainties will likely persist. 

B. Thermal photons 

The “thermal photons” are those photons result¬ 
ing from the interaction of thermalized medium con¬ 
stituents^. The computation of photon production rates 
may be done using thermal held theory techniques, or 
using relativistic kinetic theory [51]. Both approaches 
have contributed to the compendium of rates used in this 
work. 

In the partonic sector, photon-production processes 
calculated at leading order in the strong coupling con¬ 
stant, ps, have been available for almost 15 years [52]. 
Those are used here^. At high energies, the charged 
particle multiplicity is dominated by mesons. In the 
hadronic sector at temperatures comparable to, and 
lower than, the crossover temperature, photons originat¬ 
ing from thermal reactions of mesonic origins were calcu¬ 
lated in Ref. [26] . That same work also includes the pho¬ 
tons obtained from taking the /9-meson self-energy to zero 
invariant mass. This procedure accounts for the baryonic 
contributions, be it radiative decays or reactions of the 
type ttN —>■ ttN'j, and NN —>■ NN'y, where N represents 
a nucleon. The net rate parametrized in Ref. [26] also 


^ The thermalization approximation will be relaxed later. 

^ Some recent work has extended this seminal result by going up 
next-to-leading order [53]. For values of the strong coupling rel¬ 
evant to the phenomenology considered in the current work, the 
net photon rate at NLO is a modest 20% larger than that at LO. 


As the strongly-interacting fluid hadronizes, it trans¬ 
forms into hadrons which will interact. When those in¬ 
teractions cease, the momentum distributions are frozen 
and the particles free-stream out to the experimental de¬ 
tectors. The longer-lived hadrons will contribute signif¬ 
icantly to the photon signal and therefore have to be 
included. Collectively, they are dubbed “the cocktail” 
and are (for ALICE) ,f], p,uj,r]', (j); the relevant photon- 
producing decays are subtracted from the measured in¬ 
clusive signal [54], to expose a combination of thermal 
photons and prompt photons. There are however other, 
shorter-lived, states which decay with a photonic com¬ 
ponent in the final states [55]. This work includes all 
of the ones with a mass M < 1.7 GeV. The differential 
cross section of the decay photons can then be calculated, 
knowing the relevant branching ratio. After including all 
of these, together with the decays considered in Ref. [56], 
the most important channels were found to be E —> Ay, 
/i(1285) —>■ p°y, and AT* (982) —>• ATy. All contributions 
are however included, for completeness. 

IV. CORRECTING THE PHOTON EMISSION 
RATES FOR VISCOSITY 

As mentioned earlier, it is an established fact that the 
bulk dynamics of strongly interacting matter is sensitive 
to the value of shear and bulk viscosities, two of the trans¬ 
port coefficients of QCD. Switching to a corpuscular de¬ 
scription, and considering separately the reactions that, 
together, define the fluid enables a channel-by-channel 
viscous correction of the photon emission rates. The pho¬ 
ton production rate, R-y, admits a kinetic theory formu¬ 
lation. For 2 —2 scattering (1 -|- 2 —>• 3 -I- y) it is [51] 

d^Ry 1 f d^pi d?p2 d^p3 

^ d^k ~ 2(27r)3 J 2P0(27r)3 2P°{2tt)^ 

x{27t)^6\Pi+P2 - Pa - i^)|X| VB/F(i"l)/B/F(P2) 

X (l + 0 'B/f/b/f(^3)) ) (9) 

where |A4p is the squared matrix element corresponding 
to the 2 —>■ 2 scattering, /p is the particle momentum 
distribution for bosons {ap = 1) or fermions {ap = —1), 
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and the photon four-momentum is K = (w, k). The dis¬ 
tribution function must then be modified in the pres¬ 
ence of dissipative effects, for the kinetic formulation of 
to match that with the explicit dissipative trans¬ 
port coefficients. This modification is written as Jb/f = 

fs/F + ^Ib/F- 

Linearizing in S/b/f yields 


LO 


(Pk 


PrT (vise) 


OJ- 


d?k 


■ UJ- 


d?k 


( 10 ) 


where 

<PR^ 1 f Ppi Pp2 (Pps 

^ Pk ^ 2(27r)3 J 2P0(27r)3 2P°{2Tr)^ 

xi2F)^S\Pi+P2-P3-K)\M\‘^ 


X 


'SfB/FiPl)C/FiP2) (l + ^B/F CMP^)) 

+ /b/V('^i)^/s/f(-F2) (^1 + CFb/F /b/f('^3)) 
{o'b/F SfB/FiPs)) ■ 


( 11 ) 


Next, the corrections appropriate for shear and bulk vis¬ 
cosity are discussed. While there is currently a consid¬ 
erable body of research devoted to the extraction of the 
shear viscosity from relativistic heavy-ion collisions [16], 
studies of the bulk viscosity are rarer. To proceed fur¬ 
ther, assumptions need to be made about the form of 
SfB/F{Pf^)- The space-time coordinates of the emis¬ 
sion site, X, are now explicit. The first-order correction 
5fB/F{P)X) is linear in the shear stress tensor and 
the bulk pressure 11. In this case: 


5fB/F{P, X) = TT^,{X)P^P’'S{P, X) + U{X)B{P, X), 

( 12 ) 


where two properties of f^j,(X), = 0 and 

F^^{X)P = 0, were used to constrain the expansion of 
SfB/F{P,X) in TTi^^iX). 

The functions S{P, X) and B{P, X) can depend on the 
spacetime position X through e.g. the local value of the 
temperature T{X), the energy density e{X), the entropy 
density s(X), etc. All these implicit functions of X are 
thermodynamical quantities that are related through the 
equation of state of the medium. For practical reasons, it 
is better if corrections to photon emission do not have an 
explicit dependence on the equation of state. This can 
be achieved if the momentum dependence of S{P, X) and 
B{P,X) can be factorised from the rest such as: 

SfB/F{P,X) = 7r^,iX)Pf^P‘'Y.^x\x)S^M\P,T) 

j 

+ n(x)^i?«(x)i?«(p,r), (13) 

3 


where it was assumed that the momentum-dependent 
factor S/B^^{P,T) could also depend on the temper¬ 
ature, but no other thermodynamical quantities. The 


subscript X was used to identify the spatial part of S 
and i?, while the subscript M is used for the momen¬ 
tum dependent term. The sum over j is necessary if 
e.g. B{P,X) cannot be factorised as Bs{X)Bm{P) but 
is factorisable as a sum of such terms {B‘g\x)B‘'^'^ (P) + 
B'g\x)B^^\P)). It will be shown shortly that such a 
general form is appropriate for the 6fB/F{P,X) used in 
this work. More generally, this factorization and the ex¬ 
pansion in Eq. (13) can be seen as an expansion of irre¬ 
ducible tensors in momentum space [34]. 

Using Eq. (13), the effect of viscosity on photon pro¬ 
duction (Eq. (11)) can be written 

,3 „ (vise) 

U3^ =f^,{X)K>^K''Y.S%\x)S%\k,T) 

j 

+n(x)^i?«(x)B^)(if,r), (14) 

3 


where = 0 and tTi^^{X)P = 0 were used again 

to constrain the coefficient multiplying Tr^^{X). 

The coefficient S^^{K,T) is given by [57] 

3 ( 3)1 


^m{X,T) 2{K-u)2 


Qfj.i' + ‘2.u^Uv + 3 


1 


Ppi 


{K ■ uY 

Pp2 


{^K^Uu -\- UfiKY 

(K ■ u) 
d^P3 


2(2f)3 


2P0(2f)3 2P0(2f)3 2P°{2fY 
x{2fYsYPi+P2-P3-K)\M\‘^ 

X /i%(^ 2 ) (1 + CTB/Ff^B/FiP^)) 

[p^P^st}{P2)) (1 + ctb/f/b]b(P^) 

+ /B)FiPl)fB/FiP2WB/FP3^P3S^Yl\P3) 


+/b/V(^i) 


(15) 


while B^jY\k,T) is given by the simpler expression 


B^P(K T) = ^ f d^P2 d^P3 

^ > 2(2f) 3 7 2p0(27r)3 2P0(2 f)3 2P30(2f)3 

x(2f)V(Pi + P 2 - P 3 - K)\M\'^ 


X 


\b<YiHPi)) /i%(P2)(l + <JB/Fff/FiP^)) 

+-/'b/f(^i) (-®M (T 2 )) (1 + crB/F/ij/F(T3)) 

+/i%(^^i)/i°/V(^2) (aB/FB^7^3))' 


(16) 


Since S^{K,T) and B^^{K,T) are scalars, they can 
only depend on K through the combination K ■ u. It is 
thus possible to evaluate Sm and Bm in the restframe of 
the fluid where the photon energy is uj = K ■ u and the 
temperature is T. Although s')^{K,T) and b'^^{K,T) 
cannot generally be reduced to an analytical expression. 
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Rate 

Ideal 

Shear correction 

Bulk correction 

QGP — 2-^2 

[52] 

Yes [57] 

Forward scattering approximation 

QGP — Bremsstrahlung 

[52] 

No 

No 

Hadronic — Meson gas (tt, K, p, K*, ai) 

[26] 

Yes [18, 58] 

Yes [this work] 

Hadronic — p spectral function (inch baryons) 

[26, 27] 

No 

No 

Hadronic — p -\- p bremsstrahlung 

[27, 59] 

No 

No 

Hadronic — p-p-cj system 

[28] 

No 

No 


TABLE II. A summary of the thermal photon rates sources, together with the current state of advancement of their viscous 
correction. The bulk corrections are original to this work. The “forward scattering” approximation refers to considering cases 
where, in 2 —>■ 2 scattering, the exchanged momentum is soft (i.e. ~ gT ). In this case the amplitude will be dominated by 
forward scattering. 


it is nevertheless possible to tabulate them as functions 
of uj and T. 

Putting all of this together, the final expression for the 
photon emission rate is thus^ 


+ {X)Kf^K‘' Y. {K, T) 


+ n{X)YB^i\x)B^^\K,T) . (18) 

d 


Using the approximations for the hadron/parton dis¬ 
tribution functions outlined in Appendices (A) and (B), 
one may now summarize the corrections to the distribu¬ 
tion functions that arise from the inclusion of shear and 
bulk viscosity. For temperatures where the degrees of 
freedom are partonic, 


(P, X) = 7t^,{X)P^P’'Sx{X)Sm{P, T) 

+Ii{X)BY''{X)B%^^{P,T) , (19) 


^ Note that if photon production rate were to be computed using 
field-theoretical techniques, then [51] 


U) 


d^k 


1 

(27r)3 





1 ) ’ 


(17) 


where k) is the retarded, finite-temperature, photon self- 

energy. This equation is exact in the strong interaction, but 
correct to leading order in the electromagnetic coupling, a. In 
that case, linearization in can still be used to write 

the photon rate as Eq. (18), although with different expressions 
for S^P{K) and 


with (suppressing arguments) 


Sx = 


2{e + Py 


Sm = 


/b/f (l+ 


rp2 


B 


QGP 

X 


15(1 _c2)(e + P) 


B 


QGP 

M 


fn/piP) ( 


1 + 


T 

P-u 


P ■ u 


T 


( 20 ) 

For hadronic degrees of freedom (a “hadronic gas” [HG]), 
it is given by 

5f^fp{P,X) = 7r^,{X)Pf^P‘'Sx{X)SM{P,T) 

+n(A) [bY^\X)B^^’\P,T) 
+bY’\X)B^^'\P,T)\ , ( 21 ) 


with 

B 

B 

B 


HG,1 

X 

HG 


'C’ 


B 


HG,2 

X 


m _ 2 
3 " 


C 


M 

HG,2 

M 


= fs/piP) (l + ^B/pf^B/piP)'" ^ ^ 

/B/piP){^ + <^B/pf^°}BiP) 


P-u 
P ■ u 


( 22 ) 


The above decompositions are not uniquely defined, since 
temperature factors and constants can be in either coeffi¬ 
cients. This is not a problem as long as the above defini¬ 
tions are used consistently. These equations thus define 
the S and B functions that enter the calculation of the 
viscous photon rates. Then, care must be taken in eval¬ 
uating s\^\k,T) and B^^{K,T)\ a discussion relevant 
for photon production in the QGP through 2 —>■ 2 scat¬ 
tering at leading order in the strong coupling constant 
appears in Ref. ([57]). 

At present, not all photon sources known are amenable 
to a calculation of viscous (shear and bulk) corrections. 
The situation is summarized in Table II, together with 
the appropriate references. 





























V. EVALUATING THE PHOTON MOMENTUM 
ANISOTROPY 


Before comparisons with measurements are made, it is 
necessary to explain how photon momentum anisotropies 
were evaluated in this work, in view of the subtleties 
involved in their calculation. 

There is a wide literature in heavy ion physics on the 
different methods of studying the azimuthal anisotropic 
flow of final state particles (see e.g. [3] for recent review). 
The azimuthal dependence of a given particle’s underly¬ 
ing momentum distribution, in a particular event, is usu¬ 
ally characterized by a Fourier series written in the form 
of a Fourier coefficient u® and event-plane angle 4'®: 


r;®e“'^- 


f dpxdyd4>px 

'o d^N-' 
P d^p 

gin<j> 

f dpxdyd(j)px 

'o d^Np 



(23) 


The superscript “s” denotes the particle species corre¬ 
sponding to the underlying momentum distribution. It 
is possible to assign additional labels identifying the kine¬ 
matic cuts used on px and y, but this is not necessary 
in what follows, since there should not be any ambiguity 
about the cuts used for each particle species. 

Experiments measure samples of the distribution 
/(fip, averaged over numerous heavy ion col¬ 
lisions. Given an appropriate measurement of az¬ 
imuthal correlation between different particles, these 
event-averaged measurements can be mapped to event- 
averages of Vn and ff'n. These can in turn be computed 
from theoretical models which can often access the full 
p^d^N^/d^p distribution. 

Due to the limited statistics available, photon 
anisotropy measurements are not photon-photon corre¬ 
lations, but rather photon-hadron correlations. Both 
experiments that have measured photon anisotropies, 
PHENIX at RHIC and ALICE at the LHC, used the 
event-plane method [60] to make the measurement. This 
method can be understood as using hadrons to define an 
effective reference plane in the transverse direction, based 
on the hadron’s azimuthal distribution; the photon mo¬ 
mentum anisotropy is then measured with respect to this 
hadronic plane. Depending on the number of hadrons 
being measured and on the size of their azimuthal mo¬ 
mentum anisotropy, the hadronic event-plane cannot nec¬ 
essarily be reconstructed accurately. This introduces a 
small uncertainty, of the order of 10%, in the mapping 
of event-plane method measurement to the and Tn of 
photons and hadrons [61-63]. 

The two limits of event-plane method measurements 
are known as the low and high resolution limits. In the 
low resolution limit, the event-plane anisotropy reduces 
to the scalar product Vn{SP} anisotropy [3]: 


Vr^iEP} V„{SP} 


{vlv^ cos(n(ff'^ - ff-])))) 




h\2 


( 24 ) 


The other limit is the high resolution limit: 

Vr^iEP} {vl cos{n{^l - ff-O))). (25) 

The angle brackets (...) represents an average over 
events. 

The resolution correction, which quantifies the accu¬ 
racy of the event-plane reconstruction, is used to deter¬ 
mine which limit of Vn{EP'\ should be compared with 
measurements. Their value for both RHIC and LHC 
measurements are shown respectively in Refs. [12] and 
[64] . The value of the resolution correction changes with 
the centrality and methods used to determine the event- 
plane. For n = 2 {v 2 {EP}), it is neither clearly in 
the high nor low resolution limit. On the other hand, 
higher harmonics (n > 2) are closer to the low resolu¬ 
tion limit. Equation (24) is thus used in this work to 
evaluate Vn{EP}. It was verified that the other limit of 
Vn{EP}, Eq. (25), does not differ from Eq. (24) by more 
than 10% [49]. The uncertainty associated with this am¬ 
biguity in Vn{EP} is thus not a significant issue. 

The experimental measurements correlate hadrons 
from a wide bin in px to photons measured in a small 
px bin, effectively resulting in a v2{EP} differential in 
the photon transverse momentum: 

low^res. {vl{Px)Vn COs(^(^^(Pt) " ^n))) 

VM?) 

(26) 


vl{EP}{pl) 


where the h superscript refers to the charged hadrons 
with which photons are correlated. In this work we eval¬ 
uate and dtj) at midrapidity, integrated over px > 
0.3 GeV. It was verified that the result of Eq. (26) did 
not change with other choices of lower p^-cuts between 
0 and 0.5 GeV. 

Equation (26) assumes that the events that are av¬ 
eraged over have small multiplicity fluctuations. That 
is, all events are assumed to produce a similar number 
of photons and hadrons. If large multiplicity fluctua¬ 
tions are present, Eq. (26) will take a different form de¬ 
pending of the details of the measurement, for example 
whether all events are treated equally, or if events with 
more hadrons and photons are given a larger weight in 
the event-average. 

To reduce the importance of multiplicity fluctuations, 
experimental collaborations first measure v2^{EP} in 
small centrality bins [64]. The anisotropy measurements 
from these smaller centrality bins are then recombined 
into a larger centrality to reduce the statistical uncer¬ 
tainty of the measurement. When the small centrality 
bins are recombined, each centrality is weighted by the 
number of photons measured in the centrality [64]: 


u[[{EP}[Cmin,Cniax] = 


E 


CG[C: 


•min ) '-■maxj 


vZ{EP}[c]N[c] 


E 


cG[c, 


•min i'-max 


N[c 


(27) 


where A^[c] is the number of photons measured in central¬ 
ity c, v2^{EP}[c] is the momentum anisotropy measured 
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in c and [cmin, Cmax] is the final (large) centrality class 
in which the measurement is reported. At the LHC the 
sub-bins are [64] 0 - 5%, 5 - 10%, 10 - 20%, 20 - 30% 
and 30 — 40%, while 10% bins are used at RHIC [65]. 

The quantity vZ{EP}[cmin, Cmax] — Eq. (27) — is the 
one that should be compared to PHENIX and ALICE 
measurements. All photon anisotropy calculations pre¬ 
sented in this paper are computed with Eq. (27) using 
the bins just listed. 


VI. RESULTS AND DISCUSSION 

We now show and discuss the result of integrating the 
photon rates discussed in Sections III B and IV, with the 
hydrodynamic approach discussed in Section 11. Prior to 
doing this, an important clarification is needed. The 
model used here is a hybrid approach, in the sense that 
it is not purely hydrodynamics: it has a viscous fluid- 
dynamics stage that is followed by a transport phase - 
modelled with UrQMD - with dynamic decoupling. The 
UrQMD afterburner is important to a successful theoret¬ 
ical interpretation of the measured proton spectra and 
V 2 [8, 66]. However, extracting the photons via the vec¬ 
tor meson spectral density [26] from a transport model 
is still very much a topical subject of current research. 
More generally, electromagnetic emissivities are typically 
calculated in conditions near thermal equilibrium, as dis¬ 
cussed earlier in this paper, and a knowledge of the lo¬ 
cal temperature and of other thermodynamic variables is 
usually absent from most transport formulations. One 
resolution of this situation has been to coarse-grain the 
transport final states, and to assign local temperatures 
to cells on a space-time grid using the equation of state 
[67, 68]. Such procedures are numerically-intensive, but 
will be studied within our framework in detail in the fu¬ 
ture. The point of view adopted in this work is that, 
apart from proton observables, hydrodynamics does pro¬ 
vide a realistic environment for the bulk of hadronic ob¬ 
servables, especially if the bulk viscosity is included [8]. 
Therefore, for the calculation of photons, the contribu¬ 
tion of the UrQMD phase of the spatiotemporal evolu¬ 
tion is modelled by letting the fluid-dynamical evolu¬ 
tion proceed past the switching temperature from hy¬ 
dro to UrQMD (the “particlization temperature” [69]), 
Tswitch = 145 MeV, down to a more typical hydro freeze- 
out temperature of T = 105 MeV. In hydrodynamical ap¬ 
proaches in general, the freeze-out temperature is a free 
parameter of the model: more words about the depen¬ 
dence of the photon signal on this parameter will appear 
later in this section. 


A. RHIC 

The direct photon spectrum and V 2 were measured 
at RHIC by the PHENIX collaboration [10-12, 70]. 
These measurements were made in Au-Au collisions at 



Pt (GeV) 



FIG. 4. The result of a hydrodynamic calculation of direct 
photon spectra, for Au - Au collisions at RHIC, in the 0-20 
% (top panel) and 20 - 40% (bottom panel) centrality range. 
The different curves are explained in the text, the data are 
from Ref. [11]. 

^ysMN = 200 GeV for centralities 0-20% and 20-40%. 
Comparison of the hydrodynamical model’s results for 
direct photon spectra are shown in Fig. 4. The prelim¬ 
inary, minimum-bias direct photon spectrum measure¬ 
ment from STAR [71] is shown in Fig. 5, and is com¬ 
pared with both the hydrodynamical calculations and the 
PHENIX measurements from Ref. [11]. The dashed lines 
represent the thermal contributions, that is the sum of all 
contributions of thermal origin. The prompt photons are 
calculated in NLO QCD, as explained earlier. The con¬ 
tribution of non-cocktail photons (Section HIC) is also 
shown. 

The curves labeled “direct” represent the sum of all 
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Pt (GeV) 

FIG. 5. The result of a hydrodynamic calculation of direct 
photon spectra, for Au - Au collisions at RHIC, in minimum 
bias centrality range. The data are from Refs. [11, 71]. 


sources considered in this work (Section III). One ob¬ 
serves that the calculation, with the contributions enu¬ 
merated in the text, and the experimental tend to con¬ 
verge for values of px ^ 2.5 GeV. There, the calcula¬ 
tion almost entirely consists of the pQCD component. 
For intermediate transverse momenta (as defined by this 
figure, Pt ~ 1.5 GeV), the calculation underestimates 
the PHENIX data central points roughly by a factor 
of 3. Agreement of the calculations with the preliminary 
STAR data (Fig. 5) is considerably better, well within 
systematic uncertainties. 

In the low px region, calculation and data are reunited 
again, but bear in mind the strong caveats regarding the 
trustworthiness of the pQGD calculations at such low 
transverse momenta. As supported by a direct compari¬ 
son with pp photon data, the prompt photon curve shown 
in Figs. 4 and 5 should hold down to pt ~ 1 GeV. While 
one does not expect a sudden breakdown of the formalism 
used here, it does becomes less predictive as the photon 
momentum goes down. The theoretical interpretations 
of photon production in nucleus-nucleus collisions would 
rest on much firmer ground if a fundamental measure¬ 
ment of soft photons from pp collisions, extending to val¬ 
ues of transverse momenta compared to those in Figs. 4 
and 5 existed. Such a measurement, while challenging, 
would provide a valuable baseline for phenomenological 
modelling, and would further our understanding of QGD 
in its strongly coupled regime. 

Figure 6 shows the calculated photon elliptic flow, com¬ 
pared with data measured by the PHENIX collaboration. 
The photon anisotropy was evaluated with Eq. (27). The 
elliptic flow shows the now characteristic shape, with the 
turnover at pr > 2 GeV driven by the pQGD photons. 




FIG. 6. Hydrodynamic calculation of the direct photon V 2 , 
for Au - Au collisions at RHIC, in the 0 - 20 % (top panel) 
and 20 - 40% (bottom panel) centrality range . The data are 
from Ref. [12]. 


As was the case for the photon spectra the calculation of 
the photon elliptic flow systematically undershoots the 
central data points. However, and this also holds for 
the spectra, taking into account the statistical and sys¬ 
tematic uncertainties greatly reduces the tension between 
theory and experiment. Thermal photons, represented by 
the dashed curves, are shown separately to highlight that 
the thermal contribution does exhibit a large V 2 , but that 
this momentum anisotropy is then suppressed by prompt 
photons. 

As can be expected from their small contribution to 
the direct photon spectra (Fig. 4), non-cocktail photons 
do not contribute significantly to the direct V 2 ■ They are 
not shown in Figure 6. 


B. LHC 

The direct photon spectrum and V 2 in Pb-Pb colli¬ 
sions at ^/S]\fN = 2760 GeV are presented in Figs. 7 and 
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Pt (GeV) 



Pt (GeV) 

FIG. 7. The direct photon spectrum for Pb-Pb collisions at 
the LHC 0 - 20% (top panel) and 20 - 40% (bottom panel) 
centrality range. The different curves are explained in the 
text, and the data are from the ALICE Collaboration [15]. 


8 respectively. The calculations are compared with mea¬ 
surements from the ALICE collaboration [14, 15, 64]. As 
for RHIC, the contribution to the spectrum of prompt, 
thermal and non-cocktail photons, along with their sum 
(direct photons), are shown separately in Fig. 7. The el¬ 
liptic flow of thermal photons and of the total number of 
direct photons is plotted in Fig. 8. The general features of 
the photon data set at the LHC is reminiscent of that at 
RHIC, but important differences emerge when comparing 
with theoretical calculations. For both LHC observables 
— spectrum and V 2 — there is less tension between data 
and theory than at RHIC. In fact, the theory results are 



FIG. 8. The direct photon n 2 at 0 - 40% centrality. Data are 
from the ALICE Collaboration [14, 64]. 

in agreement with the experimental results when con¬ 
sidering the statistical and systematic uncertainties. As 
previously, the prompt contribution begins to take over 
at around pr ~ 3 GeV, but otherwise lies systematically 
below the thermal sources. 


C. Effect of bulk viscosity 

The calculation of direct photons presented in this 
work is the first one to include the effect of bulk viscosity 
on both the medium evolution and the photon emission 
rates. Considering that the introduction of bulk viscosity 
was shown to have a large effect on the description of the 
hadronic observables [8], it is important to highlight its 
effect on photon production. 

In Fig. 9, the direct photon spectrum and V 2 are shown 
with and without bulk viscosity for ^/snn = 2760 GeV 
Pb-Pb collisions at the LHC. Since the inclusion of 
bulk viscosity modifies the shear viscosity necessary to 
describe the hadronic momentum anisotropies [8], two 
direct photon calculations without bulk viscosity are 
shown: one with rj/s = 0.095, which is the shear viscos¬ 
ity necessary to describe the hadronic Vn in the presence 
of bulk viscosity [8], and one with rj/s = 0.16, for which 
a good description of hadronic can be achieved with 
C/s = 0. It can be seen that the two calculations that 
do not include bulk viscosity are similar, both for the 
spectrum and the V 2 - 

The effect of bulk viscosity on the spectrum of direct 
photons is small, and consists of a slight softening of the 
spectrum. Like the spectrum, the V 2 increases at low pT 
and decreases at high px- This changes the shape of U 2 , 
whose maximum value is shifted toward lower px- This 
is a distinctive photonic signal of the finite bulk viscosity 
of QCD around the transition region. 

The effect of bulk viscosity on direct photons can be 
divided in two separate contributions: its effect on the 
photon emission rates (Section IV), and its effect on the 
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FIG. 9. Effect of bulk viscosity on the direct photon spec¬ 
trum (top panel) and V 2 (bottom panel) in Pb-Pb collisions 
at y^sivTv = 2760 GeV. ALIGE measurements [13, 14, 64] are 
shown for reference. 


spacetime evolution of the medium. The effect of bulk 
viscosity on the emission rates is illustrated in Fig. 10 
by showing the photon spectrum and V 2 with and with¬ 
out corrections to the rates due to bulk viscosity. The 
effect of the shear viscous correction to the photon rates 
is shown as well, for reference. 

Viscous corrections to the rates have a small effect on 
the direct photon spectrum. This can be understood 
from the fact that viscous corrections are larger at higher 
Pt, where prompt photons dominate over thermal ones. 

The direct photon U 2 , on the other hand, is suppressed 
at higher pT by both shear and bulk viscosity corrections 
to the photon rates. The suppression is of the order of 
20 — 30%. Recall however that not all photon emission 
rates are corrected for the effect of shear and bulk vis¬ 
cosities, as listed in Table II. In consequence, the results 
shown in Fig. 10 most likely underestimate the effect of 
viscosity on the photon rates. 

The effect of bulk viscosity on the spacetime descrip¬ 
tion of the medium is illustrated in Fig. 11. The change 


FIG. 10. Effect of viscosity corrections to the photon emis¬ 
sion rates for the direct photon spectrum (top panel) and V 2 
(bottom panel) in Pb-Pb collisions at y/SNN = 2760 GeV. 


in spacetime volume induced by the inclusion of bulk vis¬ 
cosity is shown for different ranges of temperature on the 
left, while the effect of bulk viscosity o n the flow velocity 
distribution, as quantified by = -y/l -I- (u^)^ -|- (it^)^, 
is shown on the right. The effect of bulk viscosity is 
clear: it reduces the transverse expansion of the medium 
at low temperature, but considerably increases its space- 
time volume. Since thermal photon emission is propor¬ 
tional to the spacetime volume, the increase in volume 
translates into a larger number of emitted photons. On 
the other hand the slower transverse expansion implies a 
softer photon spectrum, with more soft photons emitted 
but less hard ones. It is the combination of these two 
effects that produce an overall softening of the photon 
spectrum in the presence of bulk viscosity. 


D. Effect of photon emission rates 

Calculations of direct photons in heavy ion collisions 
do not always use the same photon emission rates in the 
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FIG. 11. Event-average spacetime volume (left) 

and event-average flow velocity {u^)t (right) for hydrodynam- 
ical model with and without bulk viscosity in Pb-Pb collisions 
at y/SNN = 2760 GeV. 
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evaluation of thermal photons, which has a large im¬ 
pact on the level of agreement with data. The photon 
rates used in this work were summarized in Table II. The 
contribution to photon emission of a n-p-uj system [28] 
has just been published and was not included in previ¬ 
ous calculations of direct photons. More importantly, 
parametrizations for the photon emission rate evaluated 
with the p spectral function, along with additional emis¬ 
sion from TT -h TT bremsstrahlung, were made available in 
Ref. [27]. In consequence, calculations of direct photons 
made before this point often included only photon emis¬ 
sion from a meson gas. The importance of the different 
hadronic photon emission channels on photonic observ¬ 
ables is shown in Fig. 12. Since the effects of shear and 
bulk viscosity have not yet been evaluated for all these 
photon emission rates, corrections to the photon rates 
due to viscosities are not included for any emission rates 
in this comparison. 

It is clear from Fig. 12 that including only photon emis¬ 
sion from a gas of mesons leads to a considerable under¬ 
estimation of the direct photon V 2 - The photon channels 
evaluated with the p spectral function are especially im¬ 
portant. 

On a last note, it is relevant to highlight that QGP and 
hadronic photon emission rates that are altogether differ¬ 
ent from those used in the present work have been investi¬ 
gated over the past years. For completeness, the results 
of folding the hydrodynamical description of heavy ion 
collisions presented in this work with two of these rates 
are presented. The first rate is the “semi-QGP” photon 
emission rate [24], which includes confinement effects on 
photon emission. The second rate is the hadronic rate 
from Zahed and Dusling [72], which is evaluated using a 
different approach than the hadronic rates used in this 
work. Once again, viscous corrections to the photon rates 
are not included in this comparison. 

The semi-QGP photon rate is considerably smaller 
than the QGP rate, which results in a 30% suppression of 
the direct photon spectrum, shown on the top of Fig. 13. 
The U 2 , shown on the lower part of the figure, does not 
change significantly: an intuitive way of understanding 



FIG. 12. Importance of different hadronic photon production 
channels on the direct photon spectrum (top) and V 2 (bottom) 
in Pb-Pb collisions at ^/snn = 2760 GeV. 


this result is to note that while a suppression of the pho¬ 
ton rate at high temperature will increase the thermal 
photon V 2 , it will reduce the contribution of thermal pho¬ 
tons with respect to prompt photons. The two effects 
largely cancel out. 

An important consequence of the suppression of the 
QGP rate studied in Ref. [24] is that it does not match 
well anymore the hadronic rate in the deconfinement re¬ 
gion, as was the case with the QGP rate used previously 
in this work (Fig. 3). This fact has yet to be addressed 
in a satisfactory fashion; it highlights the importance of 
understanding the photon rates in the transition region. 

The hadronic rate from Ref. [72] is around 40%-100% 
larger than the one used in the present work. This results 
in a larger direct photon spectrum and a larger U 2 , as il¬ 
lustrated in Fig. 13. Since the hadronic rates being com¬ 
pared do not include the same photon production chan¬ 
nels, this difference is not unexpected. Further studies of 
these hadronic emission rates themselves will be required 
to establish if the two approaches can be found to agree 
for comparable production channels. 
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FIG. 13. Direct photon spectrum (top) and V 2 (bottom) eval¬ 
uated with different QGP and hadronic photon emission rates, 
in Pb-Pb collisions at ^/snn = 2760 GeV. See text for details. 


E. The importance of late stage photon emission 

The current understanding of heavy ion collisions is 
that a flow velocity anisotropy is created during the 
medium’s expansion as a result of the initial anisotropy 
in the energy deposition. This anisotropy increases with 
time, although it eventually plateaus and reverses under 
the effect of the viscosities and the decreasing pressure 
gradients. Using temperature as a proxy for time, this 
development of flow anisotropies is shown at the top of 
Fig. 14 as the x-y flow asymmetry for different temper¬ 
ature ranges. The anisotropy is evaluated with respect 
to the reaction plane, with the cc-axis aligned with the 
impact parameter of the colliding nuclei. Like hadrons, 
the momentum anisotropy of thermal photons is directly 
related to this flow velocity anisotropy. This is illus¬ 
trated by plotting on the lower part of Fig. 14 the V 2 
of thermal photons emitted in different regions of tem¬ 
perature. Lower temperatures are associated with larger 
time, which in turn is associated with larger thermal pho¬ 
ton V2- 


FIG. 14. Flow anisotropy of the medium with respect to the 
reaction plane (top) and momentum anisotropy of thermal 
photons (bottom), for different temperature ranges, in Pb-Pb 
collisions at = 2760 GeV. 


The magnitude of the V 2 of photons emitted at late 
times means that they can play a large role in increas¬ 
ing the direct photon V 2 - Their importance is illustrated 
in Fig. 15 by showing explicitly the contribution to the 
direct photon spectrum and V 2 of thermal photons emit¬ 
ted between T = 105 MeV and T = 145 MeV. This 
temperature range is chosen because T = 145 MeV is 
the switching temperature between hydrodynamics and 
UrQMD that best describes the hadronic observables [8] 
in Pb-Pb collisions at ^/snn = 2760 GeV, as explained 
at the beginning of this section. Photons emitted below 
this switching temperature would thus be best evaluated 
using UrQMD [73], but this challenging task will be ad¬ 
dressed in future work. The large contribution to the 
direct photon V 2 of these late stage photons highlights 
the importance of further studying photon emission dur¬ 
ing this phase of the evolution. 
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FIG. 15. Importance of post-particlization photon produc¬ 
tion on the photon spectrum (top) and V 2 (bottom) in Pb-Pb 
collisions at ^snn = 2760 GeV. 


VII. CONCLUSIONS 

In this paper, the production of photons in heavy ion 
collisions was studied at RHIC and the LHC using a hy- 
drodynamical model of heavy ion collisions. This com¬ 
prehensive model included realistic initial conditions (IP- 
Glasma), along with second-order hydrodynamics equa¬ 
tions with both shear and bulk viscosities. 

With the inclusion of the elements comprising this pa¬ 
per, most direct photon theoretical results were found to 
lie either within the limits set by the systematic and sta¬ 
tistical uncertainties of measurements at both colliders, 
or slightly below. A larger discrepancy between theory 
and PHENIX data remains, but in all cases the agree¬ 
ment between fluid dynamical calculations and exper¬ 
imental data were found to be improved compared to 
what they were in the past. This level of agreement with 
data provides strong support to the idea that thermal 
photons are the principal source of the low pT direct 
photon enhancement and of the large photon momentum 
anisotropy. 


The presence of bulk viscosity in the hydrodynamical 
evolution produced a small effect on the photon spec¬ 
trum at the LHC, and a modest change in the overall 
magnitude of the photon v^- On the other hand, it in¬ 
duced a clear change in the shape of V 2 , enhancing it at 
low p7p and reducing it at higher p7p. While theoretical 
and experimental uncertainties do not currently permit 
to determine if this change in the shape of V 2 is favored 
by data, the reduction of both uncertainties in the future 
could allow direct photons to be used to constrain the 
bulk viscosity of QCD. 

The significant contribution of late stage photon emis¬ 
sion to the photon V 2 , quantified in Section VIE, high¬ 
lights the need for a more sophisticated study of photon 
emission in this phase of the collisions. Work also remains 
to be done on constraining the photon emission rates, 
especially in the deconfinement region. These questions 
will be addressed in future work, and may help shed light 
on the different level of agreement observed with RHIC 
and LHC data. 

Another topic necessitating further attention is the 
need for a more sophisticated treatment of prompt pho¬ 
ton production in heavy ion collisions that includes both 
parton energy loss and jet-medium photon production 
[74]. This could be especially important at RHIC, where 
thermal photons are not as dominant over prompt pho¬ 
tons as at the LHC. The production of photons during 
the pre-thermalised phase of the collisions is also a part of 
the framework needing greater scrutiny. The encourag¬ 
ing results presented in this paper signal that the current 
understanding of thermal photons is mature enough for 
such investigations to be undertaken. 
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Appendix A: Chapman-Enskog theory in the 
relaxation time approximation 


The first-order Chapman Enskog approximation leads 
to the following formula for the leading nonequilibrium 
correction to the distribution function, (5/ik7 






E, 


tr 


(Al) 


where D = and The collision term 

was already simplified using the relaxation time approx¬ 
imation C [J/jk] = -ElSfik/TR, with El = and tr 
the relaxation time. We further defined the local equi¬ 
librium distribution function which corresponds to 
Fermi-Dirac or Bose-Einstein distribution. 

If the chemical potential and its derivatives are always 
set to zero, as assumed in this paper, the derivatives of 
the local equilibrium distribution function can be written 
as 


= -^/ok/ok (-E^V^lnr + ,(A2) 
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and (5/ik can be expressed in the following way, 
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Then, using the well known thermodynamic relations 
(at vanishing chemical potential) 
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and the conservation laws, up to first order in Knudsen 
number. 
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We then obtain the following expression for Sfik 


viscous pressure, 11 
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leading to the usual first order relations, 
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where we introduced the transport coefficients fj and C, 
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Using the first order relations derived above, = 
2riTRa^'' and 11 = —(tr9, to replace and 9 by 
and n, respectively, we obtain our final expression for Sf, 


Sftk = /ok/ok . ^ T^^k^k'' 
11 
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Thus, the bulk correction to the single particle distri¬ 
bution function becomes. 




ik 


“ji^n/ok/ok 



(A9) 


This is the expression used in the hadronic phase, both 
for photon and hadron emission, with the index “i” cor¬ 
responding to a given hadron or resonance. When calcu¬ 
lating (C, all hadrons and resonances with masses up to 
2.25 GeV were included. 


Appendix B: Chapman-Enskog theory in the 
relaxation time approximation with thermal masses 


dfi-k 

Tr 


- V K 


rp /ok/ok 




The shear and bulk viscosity coefficients can be iden- 
tihed by replacing the derived formula for J/^k into the 
definitions for the shear stress tensor, and the bulk 


Bartons carry thermal masses. This changes some of 
the steps in the derivation presented in the previous Ap¬ 
pendix. In this case, we shall only keep the terms that 
contribute to the bulk viscosity corrections to d/ik- We 
consider the effective kinetic theory of quasiparticles de¬ 
rived in Ref. [75], 
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For our purposes, it is enough to assume that the mass 
goes with the temperature as = gi (T) T. Then 


9m, 9/,k _ 5/,k 


As before, the leading order solution originating from 
the Chapman-Enskog approximation satisfies. 


+ A:fV^/ok - w. 


dnij ^ dj^ 
9x 3k 


O'" = -A, 


, <5/,k 


Since the mass depends on the temperature, the deriva¬ 
tives of the local equilibrium distribution function have 
to be re-evaluated and one obtains that (at vanishing 
chemical potential) 
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In addition. 
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Using these results, we obtain the following expression 
for Sfik 
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Finally, using 11 = —(9, the bulk 6f can be written as 
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where we used that tr = rn, with rn being the bulk 
viscous pressure relaxation time. 

For partonic degrees of freedom, the ratio is ap¬ 
proximated as 
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which was derived in the quasi-conformal limit [35]. 
Using the above, we obtain 
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If the running of gs is neglected, as is the case in the 
present work, the expression reduces to 


Sftp = -/^k/ok 


^-A,k 

n 


3 


-E,i,D\iiT + ^A^pkikf9 + 


3U 


iU. 


E,k dT 


15Tie + V){P-ciy 


(B7) 


We note that for a quark, anti-quarks, and gluons, the 
(B3) values of the asymptotic masses are 


Above, we only kept the terms which contribute to the 
bulk viscosity correction, i.e., terms that are scalars in 
momentum space. Next, using the conservation laws up 
to first order and thermodynamic relations, we can fur- 
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